setwd("C:/replication")
library(xlsx)
library(ggplot2)
library(reshape2)
library(ggpubr)
library(dplyr)
library(SDMTools)
library(ggrepel)
library(deldir)
library(psych)

###fig1
modello <- read.xlsx("models.xlsx", 
                     sheetName = "Foglio1")

modello$vars <- factor(modello$vars, levels = c("Returnability", "Pre-electoral coalition",
                                                "Ideological range (left-right)", "Number of parties",
                                                "Minimal winning coalition", "Majority government (in both chambers)"))

ggplot(modello, aes(vars, estimate)) +
  theme(axis.title.x =element_text(margin=margin(0,20,0,0)))+
  theme(axis.title.x = element_blank())+
  geom_hline(yintercept=0, lty=2, lwd=.5, colour="grey50") +
  geom_errorbar(aes(ymin=lowerCI, ymax=upperCI), 
                lwd=.7, width=0, position = position_dodge(width = 0.5)) + 
  labs(x ="")+
  geom_point(size=3, pch=20,position = position_dodge(width = 0.5)) +
  coord_flip() +
  xlab("Variable") +
  scale_y_continuous(limits = c(-9, 9), name = "Coefficient") +
  theme_minimal() +
  theme(text = element_text(size=14))+
  facet_grid(~model)

###fig2

diff <- read.xlsx("predictive performance1.xlsx", 
                  sheetName = "Foglio1")


diff$gov2 <- factor(diff$gov, levels = unique(diff$gov))
levels(diff$gov2) <- gsub("\\(", "\n(", levels(diff$gov2))

ggplot(diff, aes(x = gov2, y = diff, fill = Model)) + 
  theme_minimal() + 
  geom_bar(stat = "identity", width = 0.5, position='dodge') +
  theme(axis.text.x = element_text(angle = 65, hjust = 1)) +
  theme(text = element_text(size=14)) +
  xlab("Government") +
  ylab("Difference in predicted choice probability") +
  theme(legend.position="bottom") +
  scale_fill_grey()

###fig3

prob <- read.xlsx("predictive performance2.xlsx", 
                  sheetName = "Foglio1")


prob$gov2 <- factor(prob$gov, levels = unique(prob$gov))
levels(prob$gov2) <- gsub("\\(", "\n(", levels(prob$gov2))

ggplot(prob, aes(x = gov2, y = prob, fill = Model)) + 
  theme_minimal() + 
  geom_bar(stat = "identity", width = 0.5, position='dodge') +
  theme(axis.text.x = element_text(angle = 65, hjust = 1)) +
  theme(text = element_text(size=14)) +
  xlab("Potential government") +
  ylab("Predicted choice probability") +
  scale_fill_grey() +
  theme(legend.position="bottom")

###fig 4

esperti2018 <- read.xlsx("esperti2018.xlsx", 
                  sheetName = "Sheet1")

data_mean <- group_by(esperti2018, party, dimension, scale, vote_share) %>% summarize(mean = mean(score, na.rm=TRUE),
                                                                       sd = sd(score, na.rm=TRUE))
importance <- subset(data_mean, scale == "importance")

data_importance <- group_by(importance, dimension) %>% summarize(wmean = weighted.mean(mean, vote_share),
                                                                 sd = wt.sd(mean, vote_share),
                                                                 N= length(mean),
                                                                 se   = sd / sqrt(N))

data_importance$error <- qnorm(0.975)*data_importance$sd/sqrt(data_importance$N)
data_importance$lowerCI <- data_importance$wmean-data_importance$error
data_importance$upperCI <- data_importance$wmean+data_importance$error

ggplot(data=data_importance, aes(x=reorder(dimension, +wmean), y=wmean)) +
  geom_bar(stat="identity", width=0.5,fill=c("#999999") ) +
  coord_flip() +
  ylab("Salience") +
  xlab("Dimension") +
  theme_minimal() +
  theme(text = element_text(size=16)) +
  geom_errorbar(aes(ymin=lowerCI, ymax=upperCI),
                width=.2, position=position_dodge(.9)) +
  scale_y_continuous(limits = c(0, 20))

###table 3 (factor analysis)

position <- subset(esperti2018, scale %in% c("position"))
position_wide <- reshape2::dcast(position, id + party + vote_share ~ dimension, value.var="score")

position_wide.active <- position_wide[,c(4:10, 12, 14)]

fit <- fa(position_wide.active, nfactors = 2, fm="minres", rotate="varimax", 
          weight = position_wide$vote_share)

fl <- fit$loadings[,]
print(fl)

###fig5

scores <- data.frame(fit$scores)

names(scores) <- c("dim1", "dim2")

data_fa <- data.frame(position_wide[,c(1:3)], scores)

data_fa_long <- reshape2::melt(data_fa, id=c("id", "party", "vote_share"), measure.vars = c("dim1", "dim2"))

data_fa_map <- group_by(data_fa_long, party, variable, vote_share) %>% summarize(mean = mean(value, na.rm = TRUE))

map_wide <- reshape2::dcast(data_fa_map, party + vote_share ~ variable, value.var="mean")

voronoi <- deldir(map_wide$dim1, map_wide$dim2)

fig5 <- ggplot(map_wide, aes(x = dim1, y= dim2)) + 
  geom_text_repel(aes(label=party, size = vote_share), 
                  nudge_x = 0, nudge_y = 0, xlim = c(-1.8,1.8), ylim = c(-1.8,1.8), direction = c("y")) +
  scale_size(range=c(2,8), name = "Vote share") +
  scale_x_continuous(limits = c(-1.8, 1.8), name = "Economic left-right") +
  scale_y_continuous(limits = c(-1.8, 1.8), name = "Integration-Demarcation") +
  geom_vline(aes(xintercept=0), linetype="dashed") +
  geom_hline(aes(yintercept=0), linetype="dashed") +
  theme_minimal() +
  geom_segment(
    aes(x = x1, y = y1, xend = x2, yend = y2),
    size = .5,
    data = voronoi$dirsgs,
    linetype = 2) +
  theme(legend.position="bottom") +
  theme(text = element_text(size=16))

x <- data.frame(map_wide$dim1, map_wide$dim2)
dist <- dist(x, method = "euclidean")


coord1 <- data.frame(x1=-0.6862290, y1=0.81990751, x2=0.7827041, y2=0.90669160)
coord2 <- data.frame(x1=0.7827041, y1=0.90669160, x2=-0.3440175, y2=-0.87784494)
coord3 <- data.frame(x1=-0.6862290, y1=0.81990751, x2=-0.3440175, y2=-0.87784494)

fig5 + geom_segment(data=coord1, aes(x=x1, y=y1, xend=x2, yend=y2)) +
  geom_segment(data=coord2, aes(x=x1, y=y1, xend=x2, yend=y2)) +
  geom_segment(data=coord3, aes(x=x1, y=y1, xend=x2, yend=y2)) +
  annotate("text", x = -0.1, y = 1, label = "1.471") +
  annotate("text", x = -0.65, y = -0.2, label = "1.732") +
  annotate("text", x = 0.3, y = -0.2, label = "2.110")
